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METHODOLOGY  OF  AN  EVENT-DRIVEN  MONTE  CARLO  MISSILE  SIMULATION 


Mary  Robin  Holliday 
Center  for  Naval  Analyses 
4401  Ford  Avenue,  PO  Box  16268 
Alexandria  VA,  22302-0268 


Abstract.  This  paper  describes  the  methodology  used  by  the  Center  for  Naval  Analyses  to  develop  a 
Monte  Carlo  missile  simulation  that  computes  probabilities  of  target  acquisition  and  distributions  of 
missile  arrival  times.  The  simulation  accounts  for  numerous  missile  and  target  uncertainties. 
Coordinated  strikes  are  simulated  by  modeling  multiple  missiles  from  multiple  launch  points  firing 
on  a  group  of  targets.  Tactical  applications  of  the  simulation  required  that  it  be  computationally  effi¬ 
cient.  This  led  to  an  event-driven  time  advancement  scheme. 


INTRODUCTION  AND  OVERVIEW 

The  Center  for  Naval  Analyses  had  developed  a  Monte  Carlo 
missile  simulation  that  computes  probabilities  of  target 
acquisition  (PACQA)  and  distributions  of  missile  arrival 
times.  The  simulation  accounts  for  numerous  missile  and 
target  uncertainties.  Targets  may  have  independent  motion  or 
correlated  motion  and  are  permitted  to  vary  direction  and 
speed.  A  Random  Tour  motion  model  is  used  to  determine 
the  times  a  target  changes  direction.  Coordinated  strikes  are 
simulated  by  modeling  multiple  missiles  from  multiple  launch 
points  firing  on  a  group  of  targets.  Each  missile  may  have  a 
unique  combination  of  flight  path  and  search  mode  and  may 
be  launched  at  any  specified  time  during  the  engagement 

Tactical  applications  of  the  simulation  include  engagement 
planning  and  determination  of  salvo  size  based  on  prob¬ 
abilities  of  acquisition.  Distributions  of  missile  arrival  times 
can  aid  in  planning  coordinated  strikes  and  in  indicating  if 
target  defenses  can  be  saturated. 

The  simulation  is  event  driven.  Parameters  are  sampled 
each  Monte  Carlo  iteration,  and  the  times  of  their  occur¬ 
rences  are  determined.  Critical  events  are  those  whose 
times  of  occurrence  determine  time  steps  in  the  simulation. 
Each  such  time  is  defined  as  a  time  of  event  (TOE).  All 
TOEs  are  sorted  in  increasing  order  and  the  differences 
between  successive  TOEs  are  used  as  the  time  steps  for  a 
Monte  Carlo  iteration. 

Once  missile  and  target  parameters  are  known  (by  Monte 
Carlo  sampling)  and  the  TOEs  are  determined,  it  is  possible 
to  determine  which  target,  if  any,  would  be  acquired.  A 
transcendental  expression  relating  target  and  missile  motion 
is  used  to  determine  the  times  when  targets  move  into  the 
seeker  swath  and  their  corresponding  locations.  If  multiple 
targets  meet  the  detection  criteria,  then  the  missile  seeker  is 
modeled  to  determine  which  target  was  acquired  first.  An 
iteration  is  complete  when  all  missiles  either  acquire  a 
target  or  complete  their  search  unsuccessfully.  This  event- 
driven  time  advancement  allows  the  simulation  to  maintain 
fidelity  and  still  be  computationally  efficient. 


SIMULATION  APPLICATIONS: 

Consider  the  scenario  in  figure  1.  The  launch  point  is 
defined  as  the  geographic  location  from  which  N  missiles  are 
launched.  The  intended  target  is  the  unit  which  a  missile  is 
launched  to  acquire.  The  Area  of  Uncertainty  (AOU) 
expresses  the  uncertainty  in  target  location.  Target  positions 
are  assumed  to  be  distributed  according  to  a  bivariate  normal 
distribution,  typically  a  90  percent  containment  ellipse.  The 
AOU  may  reflect  both  errors  in  sensor  localization  and 
effects  of  target  motion.  The  objective  of  the  scenario  is  to 
conduct  a  coordinated  missile  engagement  against  the  in¬ 
tended  targets.  In  figure  1,  the  dashed  lines  represent  the 
missile  flight  trajectories.  The  “springs”  between  targets 
indicate  correlated  target  motion  where  the  lead  unit  is 
located  inside  the  hexagon. 


FIG.  1.  Schematic  of  a  coordinated  engagement 


The  simulation  computes  the  Probability  of  target  Acquisi¬ 
tion  (PACQA)  and  missile  arrival-time  distributions.  These 
distributions  are  calculated  for  a  given  missile  acquiring  a 
given  target  as  well  as  for  multiple  missiles  acquiring  a 
common  target  PACQA  is  calculated  while  considering  all 
units  near  the  target:  hostile,  neutral,  and  friendly.  The 
probability  of  acquiring  each  unit  in  the  area  is  determined 
(background  units  as  well  as  the  intended  target).  PACQA 


is  not  the  probability  the  missile  will  hit  or  kill  the  target 
Specifically,  PACQA  does  not  consider  missile  reliability, 
survivability,  or  lethality. 

The  time  each  missile  arrives  on  a  target  makes  up  one 
point  of  a  notional  arrival-time  distribution.  The  arrival- 
time  distributions  are  reported  by  means  of  bars  whose 
widths  represent  time  intervals  and  whose  height  are 
proportional  to  the  total  number  of  arrivals  associated  with 
those  times.  Therefore,  an  arrival-time  distribution  can  aid 
in  determining  the  most  likely  time  the  missile  will  arrive 
on  target.  Arrival- time  distributions  can  also  be  generated 
for  multiple  missiles  arriving  on  a  common  target  and  to 
indicate  if  it  is  possible  to  saturate  target  defenses.  For 
example,  let  AT  represent  the  time  interval  to  achieve 
saturation  of  target  defenses.  Let  K  represent  the  number 
of  missiles  that  must  arrive  within  time  interval  AT  to 
achieve  saturation.  It  is  therefore  possible  to  predict  if  this 
saturation  occurs  and  the  most  likely  time  of  occurrence. 

ASSUMPTIONS  AND  DESIGN  CONSTRAINTS 

The  simulation  was  designed  to  be  used  as  a  “real-time”  deci¬ 
sion  aid  by  the  tactical  commander  located  at  the  launch  point 
Input  parameters  are  entered  by  an  operator  and  the  simulation 
output  is  required  within  several  minutes.  Therefore,  the  sim¬ 
ulation  was  designed  to  be  computationally  efficient  In 
addition,  the  simulation  must  run  on  a  desk-top  computer 
since  all  launch  points  may  not  have  large  computer  facilities. 
Also,  the  simulation  must  capture  the  physics  of  missile  and 
target  motion  and  their  associated  uncertainty  with  as  much 
fidelity  as  possible  given  the  first  two  constraints.  The 
PACQAs  and  missile  arrival-time  distributions  calculated  by 
the  simulation  must  be  of  high  enough  quality  to  assist  a  tac¬ 
tical  decision-maker. 

The  simulation  captures  the  physics  of  the  problem  by  using 
the  actual  missile  flight  program  to  generate  a  nominal 
search  plan  (algorithm  specific  to  a  given  missile  type). 
Monte  Carlo  is  used  to  account  for  missile  and  target  uncert¬ 
ainties.  Since  targets  rarely  have  constant  direction  and 
speed.  Random  Tour  is  used  to  model  target  motion.  The 
simulation  maintains  speed  by  using  an  event-driven  meth¬ 
odology  and  computationally  efficient  sorting  and  sampling 
techniques.  In  addition,  the  solution  technique  reduces  the 
number  of  unknowns  such  that  only  a  transcendental  expres¬ 
sion  is  solved  each  Monte  Carlo  iteration.  The  simulation 
was  developed  on  a  VAX  and  is  currently  being  transferred 
to  a  desk-top  computer.  The  simulation  requires  only  a 
random-number  generator. 

SIMULATION  STRUCTURE 

The  simulation  can  be  divided  into  six  major  sections,  each  of 
which  defines  a  computational  task  (see  table  1).  The  input 
parameters  to  the  simulation  and  the  determined  nominal 
missile  behavior  are  constant  throughout  an  engagement; 
therefore,  these  “front-end”  calculations  occur  only  once.  In 
the  latter  four  sections  of  the  simulation,  the  sampling  of 
missile  and  target  parameters  and  the  subsequent  deter¬ 
mination  of  PACQA  are  performed  every  iteration.  Dis¬ 
cussions  of  each  of  the  six  sections  follow. 


TABLE  1  Simulation  structure 


Calculations  1.  Input  parameters 

performed  once  2.  Design  of  the  missile  search  plan 


Calculations  3.  Monte  Carlo  sampling  of  missile 

performed  every  and  target  parameters 

iteration  Time  4.  Definition  of  Critical  Events 

5.  Determination  of  Intercept  Time 

6.  Target  acquisition 


Input  Parameters 

The  input  parameters  of  the  simulation  are  given  in  table  2. 
These  input  parameters  describe  the  overall  engagement — 
that  is,  the  number  of  launch  points,  corresponding  salvo 
size,  and  intended  target.  More  detailed  information  is 
required  to  plan  the  flight  trajectory  a  missile  will  use  to 
search  for  its  intended  target  In  table  2,  a  target’s  radar 
cross  section  (RCS)  is  reflected  through  the  target’s  detec¬ 
tion  range.  A  target  must  be  located  in  the  missile’s  search 
area  and  seeker  swath  as  well  as  within  its  detection  range 
in  order  to  be  acquired.  The  detection  range  of  a  target  is 
also  a  function  of  environmental  conditions  and  must  be 
determined  off-line  and  input  to  the  simulation. 

TABLE  2  Simulation  input  parameters 

1.  Location  of  each  launch  point  (Lat,  Long) 

2.  Launch  point  salvo  size 

3.  Launch  point  environmental  conditions  (wind, 
temperature,  rain) 

4.  Intended  target  for  each  missile 

5.  Search  mode  for  each  missile 

6.  Flight  trajectory  for  each  missile  (initial  heading,  flyout 
waypoints) 

7.  For  each  target:  a  location  (AOU  size  and  orientation), 
direction,  speed  and  associated  uncertainty,  time  of 
report 

-  Correlated  target  motion:  information  required  for 
lead  unit  only.  Uncertainties  about  formation  are 
input  for  other  units. 

8.  Target  detection  range  (n.mi) 


Location,  direction,  speed,  and  their  associated  uncer¬ 
tainties  are  required  for  each  target.  Also,  the  time-late  for 
this  information  must  be  reported.  However,  the  source  of 
the  targeting  information  is  irrelevant.  For  example,  the 
information  can  come  from  either  a  tracking  algorithm  or 
from  a  single-position  report 


Design  of  the  Missile  Search  Plan 

The  missile  search  plan  was  determined  by  coding  a  portion 
of  the  missile’s  Command  Launch  Software  (CLS)  into  the 
simulation.  The  algorithm  coded  is  specific  to  a  given  missile 
type.  The  CLS  determines  the  missile  flight  trajectory  based 
on  input  parameters  as  well  as  the  missile  search  area  (area  in 
which  missile  seeker  will  accept  returns). 

Monte  Carlo  Sampling  of  Parameters 

The  simulation  accounts  for  both  missile  and  target 
uncertainties.  Missile  uncertainties  arc  characteristic  of  the 
missile  engineering  design  and  environmental  conditions. 
The  target  uncertainties  typically  are  scenario  dependent  and 
are  treated  as  input  parameters. 

Figure  2  illustrates  some  of  the  missile  and  target  uncertainties 
included  in  the  simulation.  The  straight  solid  line  shows  the 
planned  flight  trajectory  for  a  single  missile.  The  “fork”  in  the 
path  indicates  the  position  at  which  seeker  initiation  should 
occur,  the  area  covered  by  the  seeker  (i.e.,  the  area  between 
the  parallel  solid  lines)  indicates  the  planned  area  to  be 
searched  by  the  missile. 

The  launch  point  may  not  have  perfect  navigation,  as  illus¬ 
trated  by  the  small  circle  at  the  base  of  figure  2;  therefore, 
the  launch  point  actually  may  be  several  miles  from  the 
estimate.  This  error  could  be  important,  because  the  de¬ 
signed  missile  flight  trajectory  is  based  on  target  and  launch 
point  latitudes  and  longitudes.  The  bearing  and  range  to  the 
target  is  calculated  from  this  information.  As  is  shown  in 
the  figure,  launch  point  location  error  may  cause  the  missile 
to  search  a  different  geographic  location  relative  to  the 
target  position  than  was  intended.  A  missile  flight  trajec¬ 
tory  and  a  search  area  consistent  with  such  an  error  are 
illustrated  by  the  dashed  lines  in  the  figure. 

o  UNCERTAINTY  AREA 

•  HOST  PROBABLE  TARGET  LOCATION 

-  NOMINAL  MISSILE  TRAJECTORY 


NOMINAL  LAUNCH  POINT  LOCATION 


FIG.  2.  Monte  Carlo  sampled  missile  and  target  parameters 


Uncertainty  in  the  location  of  a  target  is  expressed  in  terms 
of  an  AOU.  Target  positions  are  assumed  to  be  distributed 
according  to  a  bivariate  normal  distribution.  The  AOU  is 
typically  the  90-percent  containment  ellipse  for  such  a 
distribution.  An  AOU  may  reflect  both  errors  in  sensor 
reporting  and  the  effects  of  target  motion. 

The  AOU  expands  with  time  to  account  for  target  motion. 
The  designed  missile  search  area  considers  target  motion  as 
well  as  missile  navigation  errors.  The  missile  search  area  is 
illustrated  in  figure  2  by  the  large  solid  elliptical  shape.  The 
planned  geographic  location  of  the  missile  search  area  is  less 
important  than  the  area  actually  searched.  For  example,  if  a 
missile  drifts  off  course  8  miles  during  flyout,  the  search  area 
will  be  shifted  8  miles.  This  is  illustrated  in  figure  2  by  the 
dashed  elliptical  shape. 

When  the  missile  search  pattern  is  designed,  the  nominal 
flight  trajectory  is  determined.  To  account  for  uncertainty  in 
missile  behavior,  missile  parameters  are  sampled  in  each 
Monte  Carlo  iteration.  Missile  heading  bias,  velocity  uncer¬ 
tainty,  uncertainty  in  the  missile  flight  trajectory  due  to  wind, 
and  azimuthal  drift  are  all  modeled  and  are  assumed  to  be 
normally  distributed. 

Target  location,  direction,  speed  and  corresponding  uncer¬ 
tainty  are  scenario  dependent  and  are  treated  as  input 
parameters.  The  AOU  is  used  to  obtain  the  target  location 
input  parameters.  The  AOU  is  an  ellipse  whose  center  is  the 
best  point  estimate  (BPE)  of  the  target  location.  The  user 
specifies  the  probability  that  the  target  is  in  the  ellipse.  The 
simulation  computes  the  associated  standard  deviations,  as¬ 
suming  the  target  location  is  characterized  by  a  bivariate  nor¬ 
mal  distribution.  If  independent  target  motion  is  assumed, 
then  each  target’s  location  at  the  beginning  of  an  iteration  is 
determined  by  sampling  from  this  distribution. 

The  simulation  also  accounts  for  correlated  target  motion. 
The  initial  location  of  the  lead  unit  (LU)  is  sampled  from  its 
input  AOU.  Other  units  are  located  at  their  reported  distance 
and  bearing  from  the  LU  with  some  uncertainty  to  allow  for 
station  keeping.  The  variation  in  each  unit’s  position  relative 
to  the  LU  is  an  input  parameter.  Direction  and  speed  are 
determined  for  the  LU  only.  All  other  units  are  moved 
relative  to  the  LU,  with  some  specified  uncertainty  to  further 
allow  for  station  keeping. 

The  simulation  uses  a  Random  Tour  motion  model  to  deter¬ 
mine  the  times  a  target  changes  direction.  Let  (XT)  be  de¬ 
fined  as  the  expected  number  of  direction  changes  made  by 
target  (i)  in  one  hour.  (XT)  is  based  on  observed  data  for  a 
given  target  type.  Times  between  direction  changes  are 
independent  and  exponentially  distributed  with  the  mean 
time  between  direction  changes  equal  to  (lAi)  of  an  hour. 
Target  motion  is  linear  between  direction  changes.  Ex¬ 
ponentially  distributed  random  variables  are  generated  by: 

Sk  =~h  [M1"**)]  * 


where 

S ^  =  times  between  direction  changes 
Xfc  =  uniformly  distributed  random  number. 

The  actual  direction  and  speed  of  a  target  are  determined 
based  on  the  quality  of  targeting  information,  and  the 
operator's  perception  of  future  target  motion.  All  sampling 
of  target  parameters  is  done  in  a  target’s  own  coordinate 
system,  where  the  most  probable  target  location  is  defined 
as  the  origin.  When  determining  if  a  target  lies  within  a 
missile’s  field  of  view,  the  target  coordinates  are  translated 
to  the  missile  coordinate  system  for  that  calculation. 

Definition  of  Critical  Events 

This  simulation  is  event-driven.  Parameters  are  sampled  in 
each  Monte  Carlo  iteration,  and  the  corresponding  times  of 
their  occurence  are  determined  Critical  events  are  those  whose 
times  of  occurrence  determine  time  steps  in  the  simulation. 
Each  such  time  is  defined  as  a  time  of  event  (TOE). 

A  TOE  occurs  whenever  a  missile  or  target  deviates  from  its 
previous  mode  of  operation.  The  only  TOE  associated  with 
a  target  is  the  time  of  a  direction  change.  Several  missile- 
related  events  generate  TOEs.  Missile  TOEs  include  the 
time  the  missile  seeker  turns  on,  the  time  the  missile  begins 
and  completes  searching  one  flight  leg,  and  the  time  the 
missile  seeker  turns  off.  The  missile  flight  trajectory  is 
assumed  to  be  linear  between  its  TOEs. 

Ail  TOEs  are  sorted  in  increasing  order,  and  the  differences 
between  successive  TOEs  are  used  as  the  time  steps  for  a 
Monte  Carlo  iteration.  This  event-driven  time  advancement 
allows  the  simulation  to  maintain  fidelity  and  still  be  com¬ 
putationally  efficient  An  iteration  is  complete  when  all  mis¬ 
siles  either  acquire  a  target  or  complete  their  search  unsuc¬ 
cessfully.  When  using  the  maximum  number  of  missiles  and 
targets,  as  many  as  800  times  must  be  sorted  per  Monte  Carlo 
iteration;  therefore,  a  computationally  efficient  sorting  routine 
was  required  The  shell-sort  algorithm  described  later  is  used 
in  the  simulation. 

Determination  of  Intercept  Time 

Both  missile  and  target  motion  are  linear  between  TOEs.  The 
range  at  which  a  given  target  can  be  detected  is  a  user  input 
This  range  is  limited  by  the  maximum  detection  range  of  the 
missile  seeker  and  is  a  function  of  a  target’s  RCS.  The  missile 
and  target  motion  and  the  detection  range  can  be  related  in 
one  expression  in  which  the  only  unknown  is  the  time  of 
intercept,  that  is,  the  time  at  which  the  target  first  moves  into 
the  missile’s  field  of  view. 

Each  target’s  position  and  heading  are  translated  in  eagh 
time  step  from  the  target  coordinate  system  to  each 
missile’s  coordinate  system.  The  target  coordinate  system 
is  defined  with  the  origin  located  at  the  most  probable 
target  location,  and  the  axes  are  aligned  with  the  AOU 
orientation  as  illustrated  in  figure  3.  YTAR  and  XTAR 


define  the  target  centered  coordinate  system.  The  angle  (P) 
is  defined  as  the  orientation  form  north  of  the  target- 
centered  coordinate  system. 


FIG.  3.  Target  and  missile  coordinate  systems 


YMIS  and  XMIS  define  the  missile  coordinate  system 
whose  origin  is  located  at  the  launch  point.  A  direct  line 
from  the  launch  point  to  the  most  probable  target  location 
initially  defines  the  missile  coordinate  system  as  given  in 
figure  3.  The  angle  (y)  is  defined  as  the  orientation  from 
north  of  the  missile-centered  coordinate  system.  The 
angular  difference  between  the  missile  and  target  co¬ 
ordinate  systems  is  defined  as  (a).  The  definitions  of 
additional  terms  are  given  in  table  3.  The  remainder  of 
this  section  develops  the  expressions  used  to  determine 
time  of  intercept.  The  objective  is  to  find  the  time  when 
the  missile  and  target  first  close  to  the  specified  detection 
range  L,  as  illustrated  in  figure  4. 

A  target’s  location  can  be  expressed  in  a  missile’s  coor¬ 
dinate  system  by  the  expression: 

=  T  cos  (a)+X  sin(a)  +  l£  , 

XT  =-y  sin  (a)+X  cos  (a)+X^  , 

Target  and  missile  motion  as  a  function  of  time  can  be 
expressed  as: 

xT  (t)=xTo+tVTFx{er) 

Yt  (f )  =  YTo  +tVTFy(eT) 

XM  {t)-XMo  +tVMFx  (SM) 

YM(t)  =  YMo+tVMFY(eM)  , 

As  illustrated  in  figure  4,  the  projective  distances  between 
missile  and  target  can  be  expressed  by  the  following 
equations: 

,cos¥ ,  MW)=sin¥. 


TABLE  3  Variable  Definition  for  Determination  of 
Intercept  Time 


x,r 

xt,yt 

xk‘Yk 


XM '  Y M 


L 


*To>YTo 


XMo  ’  YMo 


t 

t' 


Fx(0),Fy(0) 


x  and  y  coordinates  of  the  target  in  the 
targets  entered  coordinate  system 

x  and  y  coordinates  of  the  target  in  the 
missile-centered  coordinate  system 

projection  of  the  target  coordinates  in 
the  target-centered  coordinate  system 
onto  the  x  and  y  axis,  respectively,  of 
the  missile-centered  coordinate  system 

x  and  y  coordinates  of  the  missile  in  the 
missile  coordinate  system  (n.mi.) 

velocity  of  target  (kt) 

velocity  of  missile  (kt) 

heading  of  target  in  missile  coordinate 
system  (rad) 

heading  of  missile  in  missile  coordinate 
system  (rad) 

detection  range  of  missile  for  a  given 
target  (n.mi.) 

x  and  y  coordinates  of  the  target  in  the 
missile  coordinate  system  at  the 
beginning  of  the  time  interval 

x  and  y  coordinates  of  the  missile  in 
the  missile  coordinate  system  at  the 
beginning  of  the  time  interval 

time  (hours) 

intercept  time,  time  when  missile  and 
target  are  distance  L  apart  (hours) 

expressions  that  determine  the  correct 
trigonometric  function  to  model  target 
motion  where  0  =  0p  or  missile  motion 
where  9  =  . 


\ms 

* 


Where  cl.  -  angular  position 

of  seeker  boresight 
off  centerline 


For  each  quadrant,  there  is  an  angle  for  which  the  above 
expressions  are  satisifed.  Let  ('  be  the  intercept  time.  Sub¬ 
stitute  in  the  expressions  for  target  and  missile  motion 
where  r  —t\ 


L  cosy-XTo  +t'VTFx  (0r)-XjWc)  -t'VMFx  (@Af) 

L  siny=YTo  +  t'VTFy(®r) -  YMo  -t'VMFY  (0^)  • 

Define  the  differences  in  initial  missile  and  target  position: 

AX-Xro-XMo 
AY=rTo~YMo  . 


Define  differences  in  relative  velocity  as 


Vx  =  VtFx(S t)  -  VmFx(Qm) 
VY=VTFY(,&T)-VMFY(eM )  . 

Therefore,  the  expression  for  relative  motion  can  be 
simplified  to 

L  cos  \jf= AX  +  t'Vx 
L  sin  \|/=  AT  +  rTy  . 

Square  both  sides  of  the  equations  and  add  the  expressions 
together  to  obtain  one  expression  for  intercept  time: 


L 2  =  L2  (cos 2  y  +  sin2  y)  =  AX  2  +  AY  2 
+  t'(2AXVx  +  2AYVy)  + 1 '2  (Vx2  +  Vy2^ 

This  is  a  quadratic  in  t'  with  all  coefficients  known.  The 
quadratic  expression  will  have  no  positive,  real  solution  if 
the  target  is  not  located  in  the  area  searched  by  the  missile. 
Also,  the  intercept  time  t'  must  fall  within  the  current  time 
interval  (i.e.,  within  two  successive  TOEs)  to  be  valid. 

Target  Acquisition 

A  target  is  assumed  to  be  detected  if  it  lies  within  its 
corresponding  detection  range  and  is  located  in  the  seeker 
swath  and  missile  search  area  at  that  corresponding  time.  In 
every  time  interval,  t'  is  determined  for  each  target  and 
missile  combination.  Intercept  time  is  determined  for 
missiles  only  after  they  have  begun  their  search.  For  a  given 
missile,  the  smallest  intercept  time  determined  corresponds 
to  the  target  that  is  intercepted  first  by  that  missile.  The 
target  intercepted  first  is  not  necessarily  the  target  acquired, 
since  the  target  must  also  be  located  in  the  seeker  swath  and 
missile  search  area  at  this  time.  For  multiple  targets  whose 
intercept  time  falls  within  the  time  interval  required  for  the 
missile  to  lock  on  target,  the  missile  seeker  is  modeled  in 
more  detail. 


Once  a  missile  detects  a  target,  time  of  arrival  is  determined 
FIG  4.  Target-missile  relative  positions  assuming  a  straight-line  trajectory  from  missile  to  target 


The  target  is  assumed  to  continue  moving  after  multiple 
missile  hits  (no  mobility  kill).  An  iteration  is  complete 
when  all  missiles  either  acquire  a  target  or  complete  their 
search  unsuccessfully. 

SIMULATION  COMPUTATIONAL  EFFICIENCY 

The  simulation  is  a  temporally  event-driven  model  to 
achieve  computational  efficiency.  In  addition,  options  for 
several  underlying  processes  in  the  simulation  were  exam¬ 
ined  for  efficiency  and  accuracy  and  have  been  documented 
in  [1]  and  [2]. 

Several  techniques  for  generating  normally  distributed 
random  numbers  were  examined  for  use  in  the  simulation 

[1].  The  polar  technique,  developed  by  Atkinson  and 
Pearce,  was  chosen  [3].  The  polar  technique  is  an  improve¬ 
ment  over  the  earlier-developed  Box  and  Muller  technique 
that  generates  normally  distributed  values  given  uniformly 
distributed  values.  Box  and  Muller  maintained  a  one-to-one 
correspondence  between  the  uniformly  and  normally 
distributed  value  but  required  the  calculation  of  sine  and 
cosine.  The  polar  technique  eliminates  these  trigonometric 
calculations,  thereby  increasing  the  computational 
efficiency  of  the  algorithm.  The  polar  technique  uses  two 
uniform  random  values  to  generate  two  normal  random 
values;  however,  the  technique  rejects  values  with  a 
probability  of  (1  -  ni4)  [3].  Even  with  the  rejection  of 
variables,  the  polar  technique  was  found  to  be  between  9 
and  31  percent  faster  (depending  on  the  machine)  in  the 
FORTRAN  programming  language  than  the  Box  and 
Muller  technique  [3],  and  at  least  twice  as  fast  as 
techniques  examined  in  reference  [1]. 

The  shell-sort  algorithm  was  chosen  for  use  in  the 
simulation  [2].  The  shell-sort  required  about  4.7  CPU 
seconds  to  sort  20,000  uniformly  distributed  random 
values,  which  was  more  than  an  order  of  magnitude  faster 
than  the  examined  algorithms  [2].  The  shell-sort  algorithm 
makes  use  of  an  interchange  sort  [4],  In  the  simplest 
application  of  an  interchange  sort,  pairs  of  numbers  are 
exchanged  one  step  at  a  time  until  an  array  is  sorted.  It  is 
important  to  note  that  the  computational  efficiency  of  an 
interchange  sort  depends  entirely  on  the  array  values 
exchanged.  The  shell-sorts  enhancement  to  the  interchange 
sort  is  to  allow  array  values  to  move  greater  distances  at 


first  and  then  move  smaller  distances  as  their  final 
destination  is  approached.  The  shell-sort  is  also  referred  to 
as  the  diminishing  increment  sort.  Conceptually,  the  shell- 
sort  is  simple;  however,  the  mathematical  analysis  to 
determine  initial  and  diminishing  array  value  movement 
(jump  size)  is  complex.  This  work  uses  a  straightforward 
application  of  the  shell-sort.  The  initial  jump  size  is  chosen 
to  be  half  the  array  length.  The  jump  size  is  then  repeatedly 
halved  until  the  array  is  sorted.  The  shell- sort  algorithm  is 
programmed  in  about  20  lines  of  FORTRAN.  This 
algorithm  has  the  additional  advantage  that  only  the  array 
itself  and  one  other  value  are  required  to  be  stored. 

The  computational  efficiency  of  this  event-driven  Monte 
Carlo  missile  simulation  is  best  illustrated  by  example. 
Consider  the  scenario  in  which  two  launch  points  each  fire 
one  missile  at  their  intended  targets.  There  are  a  total  of 
six  independently  moving  targets  being  modeled.  The 
simulation,  run  for  1,000  Monte  Carlo  iterations  on  a 
VAX  1 1/8650  (with  I/O),  required  24  CPU  seconds. 
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